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Abstract. The classical single-band uncertainty model introduced by 
Bertsimas and Sim [S] has represented a breakthrough in the devel- 
opment of tractable robust counterparts of Linear Programs. However, 
adopting a single deviation band may be too limitative in practice: in 
many real-world problems, observed deviations indeed present asymmet- 
ric distributions over asymmetric ranges, so that getting a higher mod- 
eling resolution by partitioning the band into multiple sub-bands is ad- 
visable. The critical aim of our work is to close the knowledge gap on 
the adoption of multi-band uncertainty in Robust Optimization: a gen- 
eral definition and intensive theoretical study of a multi-band model are 
actually still missing. Our new developments have been also strongly in- 
spired and encouraged by our industrial partners, interested in getting 
a better modeling of arbitrary shaped distributions, built on historical 
data about the uncertainty affecting the considered real-world problems. 

Keywords: Robust Optimization, Multi-band Uncertainty, Compact 
Robust Counterpart, Cutting Planes, Probabilistic Bound. 



1 Introduction 

A fundamental assumption in classical optimization is that all data are exact. 
However, many real- world problems involve data that are uncertain or not known 
with precision. This can be due, among other things, to erroneous measurements 
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or adoptions of approximated numerical representations. It is well-known that 
neglecting such uncertainty may have very bad effects: solutions that are feasible 
may reveal to be infcasible, as a result of small deviations from nominal data, 
while solutions considered optimal may turn out to be of very low quality. Such 
risks cannot be taken in critical applications, such as contaminant detection 
in water distribution networks or the design of power grids [21]. Since the 
groundbreaking investigations by Dantzig |18j . many works have thus tried to 
find effective ways to deal with uncertainty (see [HI for an overview). 

In recent years. Robust Optimization (RO) has become a valid methodology 
to tackle uncertainty in optimization and has increasingly attracted the attention 
of practitioners, that have started many research partnerships on the topic with 
academic institutions (e.g., on optical network design [7], telecommunications 
network design [53], surgery scheduhng [IS]). The central feature of RO is to 
approach coefficient uncertainties through hard constraints, that restrict the 
feasible set maintaining only robust solutions, i.e. solutions protected against 
deviations. For an exhaustive introduction to theory and applications of RO, we 
refer the reader to the book by Bcn-Tal et al. [1] and to the recent survey by 
Bertsimas ct al. [5]. 

One of the most successful and widely- adopted RO approach is the so-called 
r'-scenario set (BS), introduced by Bertsimas and Sim in [5]. The uncertainty 
model of BS assumes that each uncertain coefficient a is a symmetric and 
bounded random variable, taking value in the symmetric interval [a — d,a + d], 
where d is the nominal value of a and d is the maximum deviation from d. 
Furthermore, the uncertainty set adopts a parameter F > to represent the 
maximum number of coefficients that are allowed to deviate from their nominal 
value. These modeling assumptions lead to a compact RO model that preserves 
the features of the original nominal model. However, the use of a single deviation 
band may greatly limit the power of modeling uncertainty, as noted, for exam- 
ple, even by Sim and colleagues in |16| . This is particularly evident in real- world 
problems, where it is common to have asymmetric probability distributions of 
the deviations, that are additionally defined over non-symmetric intervals. In 
such cases, neglecting the inner-band behavior and just considering the extreme 
values like in BS may lead to a rough estimate of the deviations and thus to un- 
realistic uncertainty set, which either overestimate or underestimate the overall 
deviation. Having a higher modeling resolution would therefore be very desir- 
able. This can be accomplished by breaking the single band into multiple and 
narrower bands, each with its own F value. This observation was first captured 
by Bienstock and taken into account to develop a RO framework for the special 
case of Portfolio Optimization [5]. It was then extended to Wireless Network 
Design [10ll7j . Yet, no definition and intensive theoretical study of a more gen- 
eral multi-band model applicable in other contexts have been accomplished. The 
main goal of this paper is to close such a gap. 

We remark that investigating the adoption of a multi-band uncertainty set 
and studying the theoretical properties of the resulting model have been dis- 
cussed with and strongly encouraged by our industrial partners, such as British 
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Telecom Italia (BT) and Nokia Siemens Networks (NSN) |28| . in past and 
present collaborations about real network design |10I11I17I23] . NSN, in partic- 
ular, has been interested in finding refined models for taking into account the 
arbitrary and non-symmetrical distributions of traffic uncertainty characterizing 
nation-wide optical networks [IHl- Our modeling and theoretical developments 
are thus also strongly based on practical industrial needs. A better modeling of 
the traffic uncertainty affecting telecommunication networks was one of the es- 
sential goals of our research activities in the German industrial research project 
ROBUKOM |11I23| . This project aimed at developing new models and algo- 
rithms for the design of robust and survivable networks, in collaboration with 
NSN and the German National Research and Education Network (DFN) |19| . 

Outline of our contributions. In this work, we study the robust counterpart 
of a Mixed-Integer Linear Program (MILP) with uncertain coefficient matrix, 
when a multi-band uncertainty set is considered. As previously noted, our main 
aim is to fill the knowledge gap on multi-band uncertainty in RO by presenting 
new theoretical results on this topic. Specifically, we: 

1. introduce a general definition of multi-band uncertainty (Section II. 1(1 : 

2. characterize a family of dominating uncertainty scenarios (Section [TTT]); 

3. derive a compact formulation for the robust counterpart of a MILP, whose 
size is smaller than that presented by us in |12I13) (due to the new domi- 
nance results) (Section [2]); 

4. define an efficient method for the separation of robustness cuts (i.e., cuts 
that impose robustness), based on solving a min-cost flow instance. In par- 
ticular, we propose an alternative proof to the one that we presented in |12| : 
the new proof takes advantage of the new dominance results and highlights 
the correspondence between the dominating uncertainty scenarios and the 
integral flows of an auxiliary min-cost flow instance (Section [3]); 

5. study the properties of a special subfamily of uncertain MILPs, namely pure 
Binary Programs with uncertainty just affecting the objective function, and 
we show that if the nominal problem is polynomially solvable, then also its 
counterpart remains polynomially solvable (Section 

6. derive a data-driven bound on the probability that an optimal robust solu- 
tion is infeasible, due to deviations that are not considered by the uncertainty 
set. We define the bound data-driven, as it exploits samples of the unknown 
arbitrary distribution of an uncertain coefficient (Section [SJ; 

We remark that our original results are not obtained by simply extending the 
proofs of [S] for single-band uncertainty, but need alternative proof strategies. For 
example, the separation of cuts imposing robustness does not anymore reduce 
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to a trivial ordering problem (see |21j). but, as we have proved, corresponds to 
solve a min-cost flow problem (see Section |2] for details). 

Moreover, our bound on the probability of robust infeasibility differs from 
those proposed in the two cornerstone RO papers [5] and [8]: in contrast to the 
bounds of these two works, that use a priori distributional information (symmet- 
ric distribution over symmetric single range), we derive a data-driven bound by 
exploiting samples of the unknown distribution. This type of bound is particu- 
larly suitable in real-world problems, where historical data about the deviations 
are commonly available. Moreover, we think that a data-driven bound is also 
more in accordance with the spirit of multi-band uncertainty: multi-band is in- 
deed particularly suitable to build strongly data-driven uncertainty sets, that 
can effectively approximate the shape of histograms built on historical data [S] . 

Finally, we note that we first presented the results about the compact ro- 
bust counterpart and the efficient separation of robustness cuts in |12j . where 
we focused attention on Linear Programs. A revised version of the work with 
the alternative proof for the separation of robustness cut is presented in |13j . 
Additionally, we remark that in this paper we have decided to focus just on 
the theoretical results. However, we have already got promising computational 
results on realistic instances of real- world problems (specifically, network design 
problems), where the uncertainty set has been defined in collaboration with our 
industrial partners. We refer the reader to [12] and [IS] for a description of these 
results. 

Review of related literature and comparisons with our work. While 
writing this paper, we became aware that a model similar to ours was presented 
in a very recent technical report (|26j). However, we remark that: 1) in our for- 
mer papers [12j . |13j and [14j . we already introduced the multi-band model and 
formalized theoretical results similar to those presented in [2B]; 2) the uncer- 
tainty model of j26| is a special case of our multi-band model, so our results are 
more general. Moreover, in contrast to our model, the model of [55] cannot shape 
any arbitrary distribution, because of its limiting assumptions (symmetric devi- 
ation range, symmetric distribution of deviations, probability of deviation that 
decreases as the value of the deviation increases). Such hypotheses, in particular 
those on symmetry, greatly reduce the ability of modeling real uncertainty dis- 
tributions and can result in excessive conservatism, as clearly pointed out even 
by Sim and colleagues in |16| . 

Reducing the conservatism of robust solutions associated with first classic 
works like [5] and [8] has been an important question approached in several 
works. We recall here some more significant references, referring the reader to 
|4I6) for a broader view. In jl6j . the conservatism is mitigated by defining new 
deviations measures, aimed at capturing asymmetrical distributions. In [5], a 
"soft robust" model is introduced by using the concept of risk measure, reducing 
the conservatism at the price of lower probabilistic guarantees of robustness. The 
concept of risk measure is also used in [3] to develop rules for constructing good 
uncertainty sets that reduce conservatism. Finally, another interesting approach 
based on the BS uncertainty model is constituted by "light robustness" 
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this approach adopts a hard upper bound on the objective value to control 
conservatism and aims at finding the most robust solution that satisfies the 
bound. 

All these works successfully reach the objective of guaranteeing a reduction 
in conservatism, proposing new interesting approaches. However, at the same 
time, we think that they miss the elegant simplicity and easiness of access that 
characterize BS and that have been undoubtedly decisive for its success also 
among practitioners and scholars from outside the Operations Research commu- 
nity. Concerning this point, we note that our work differs from the previously 
cited ones: we defined and studied a general multi-band uncertainty model that 
overcomes remarkable limitations of BS, while trying to maintain the simplicity 
of BS and offering straightforward ways to build uncertainty sets strongly based 
on available historical data. We remark that the simplicity of access and the 
data-driven characteristics have been highly requested by our industrial partners 
(e.g., Nokia Siemens Networks) in several industrial collaborations. Our model- 
ing and theoretical developments are thus also inspired by and satisfy practical 
industrial needs [25] . 



1.1 Model and Notation. 



In this work, we study the robust counterpart of Mixed-Integer Linear Program- 
ming Problems whose coefficient matrix is subject to uncertainty and the uncer- 
tainty set is modeled through multiple deviation bands. We start by considering 
a generic deterministic problem of the form: 

max ^ Cj Xj (MILP) 
Gij Xj < bi i G I 

Xj >0 j e J 

Xj e Z+ j G Jz^ J 



where / = {1,...,™} and J = {l,...,n} denote the set of constraint and 
variable indices, respectively. Assuming that the uncertainty affects only the 
coefficient matrix does not limit the generality of our study: indeed, 1) if the 
r.h.s. vector b is uncertain, we introduce an additional variable Xn+i : Xn+i = 1 
and include 6 as a column of the matrix A] 2) if the cost vector c is uncertain, 
we introduce an additional variable L, add the constraint c'x > L and change 
the objective function into maxL. 

We assume that the value of each coefficient is uncertain and is equal 
to the summation of a nominal value aij and a deviation lying in the range 
[d^ , dfj^], where d^ , df^^ e R represent the maximum negative and positive 
deviations from aij, respectively. The actual value aij thus lies in the interval 
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Wc derive a generalization of the Bertsimas-Sim uncertainty model by parti- 
tioning the single deviation band [d^j~,d^~^] of each coefficient into K bands, 
defined on the basis of K deviation values: 

-oo < < • • • < < < d° = < 4 < 4 < ■ • ■ < df/ < +O0. 

Through these values, we define: 1) a set of positive deviation bands, such that 
each band fc G {1, . . . , K~^} corresponds to the range (4~^i 4']' ^ 
ative deviation bands, such that each band k G {K~ + 1, . . . , — 1, 0} corresponds 
to the range (4~^'4] "^^"^ band k = corresponds to the single value 4 
(the interval of each band but k = is thus open on the left). With a slight 
abuse of notation, in what follows we indicate a generic deviation band through 
the index fc, with k G K ~ {K^ , . . . , —1, 0, 1, ... , K^} and the corresponding 
range by (4-1,4.]. 

Additionally, for each band k G K, we define a lower bound 1^ and an upper 
bound Uk on the number of deviations that may fall in k, with l^, Ufc G Z satisfy- 
ing < Ik < Uk < n. In the case of band fc = 0, we assume that uq = n, i.e. we do 
not limit the number of coefficients that take their nominal value. Furthermore, 
we assume that X^fes-ff ''^ — that there always exists a feasible realization 

of the coefficient matrix. 

We remark that, in order to avoid an overload of the notation, we assume that 
the number of bands K and the bounds lk,Uk are the same for each constraint 
i G I. Anyway, it is straightforward to modify all presented results to consider 
different values of those parameters for each constraint. 

All the elements introduced above define a multi-band scenario set Sm. A 
scenario S specifies the deviation dfj of each coefficient of the problem (thus, 
dij = o-ij + dfj). We say that a scenario S is feasible and belongs to the set Sm 
if and only if the deviations dfj satisfy the following three properties: 

4g[4^4^] *e/,jGJ, (1) 

h < \{jeJ: 4 G (4"'>4]}| I G /, fc G /^\{A-}, (2) 

Ik- < J:4=4-}| < uj,^ tel. (3) 

Property ([T]) enforces that the deviation of each coefficient lies in the correspond- 
ing overall deviation range. Properties (H])-® guarantee that, for each constraint 
and deviation band, the number of deviations falling in the band must respect 
the corresponding bounds. 

Before deriving a compact formulation for the robust counterpart of (MILP), 
we introduce some structural properties of the multi-band scenario set Sm ■ To 
this end, we first define the band partition of a scenario S G Sm as the family of 
subsets: 

Jifc(^) = {jG J:4G(d^\4]}, ieI,keK. (4) 

Each subset JikiS) includes the indices j G J of constraint i £ I whose deviation 
falls in band fc. It is easy to verify that such subsets define a partition of J for 
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each constraint i (i.e., Ufce/f Jik{S) — J,'^i & I, and Jik{S)riJik' {S) = 9.yi £ I, 
k,k' £ K : k ^ k'). Furthermore, we call the profile of a scenario S G Sm the 
numbers p and 9k, k Cz K defined as follows: 

fe if+ 

p = min{fc G K : li + < n} (5) 

i=K- i=k+l 



(6) 



Note that p > since uq = n and X^feeK ^fe = 

Finally, we define a dominance relation among scenarios: let S, S" G Sm be 
two feasible scenarios. Then S dominates S' if df, > df, for all i E I, j £ J. 
Obviously, a solution x that is feasible under S (i.e., '^j^ / {(iij+dfj) Xj < bi,yi G 
I) remains feasible under S' . In the following lemma, we prove that we can define 
the robust counterpart of (MILP) by limiting our attention to those scenarios 
that satisfy the following properties: 1) each deviation attains the maximum 
value of the band in which it falls; 2) the cardinality of each subset of the band 
partition ^ equals the value associated with the profile ([6]) of the considered 
feasible scenario. 

Lemma 1. Let S G Sm be a feasible scenario. 

1. If dfj ^ d^j for any j G Jik{S), i £ I , k G K, then there exists a feasible 

scenario S' £ Sm dominating S . 
^- If \ Jik{S)\ 7^ 9k for any i G I , k G K, then there exists a feasible scenario 

S' G S]\i dominating S. 

Proof. 1) Since x > and the values falling in a band are bounded from above 
by d'^j, a scenario S' such that df^ = df^ for each j G Jik{S), i G I, k G K 
dominates the scenario S. 

2) Let k^ be the highest band k G K such that \ Jik{S)\ < 9k. Then there 
exists a band fc^ < k^ such that |Jjj,2(S')| > 9^.2. Let fj, G Jik^. We define the 
feasible scenario S' by setting df^ = d^^ and — dfj for j ^ ji. Then S' 
dominates S. □ 

We call a scenario S G Sm that satisfies the properties of Lemma [T] profile and 
bound valid. Due to Lemma[Tl in the rest of the paper we can limit our attention 
to these scenarios. We now proceed to study the robust counterpart of (MILP) 
under multi-band uncertainty. 

2 A Compact Robust MILP Counterpart 

The robust counterpart of the problem (MILP) under a multi-band scenario set 
Sm can be equivalently written as: 

max Cj Xj 
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^ Xj + DEVi{x,SM) <bi i € I 

Xj >0 j £ J 

Xj G j S Jz C J, 



where DEVi(x, Sm) is the maximum total deviation induced by the scenario set 
5a/ for a feasible solution x when constraint i is considered. Due to Lemma [1] 
DEVi(a;, iSm) is equal to the optimal value of the following pure 0-1 Linear 
Program (note that in this case index i is fixed); 

DEV, {x, Sm ) = max ^ ^ x^ yfj (DEVOl) 

jGJ keK 

Y.y'i = ^eif (7) 

2/ye{0,i} jeJ,keK. (9) 

The generic binary variable (O is equal to 1 if coefficient a^j deviates in band 
k and is otherwise. The constraints ([7]) impose the profile of Sm on the num- 
ber of coefficients deviating in each band k. The constraints ([5]) ensure that 
each coefficient deviates in at most one band (actually these should be equality 
constraints, but, for the assumption uq = n made in Section ll.ll we can con- 
sider inequalities). An optimal solution of (DEVOl) thus defines a distribution of 
the coefficients among the bands which maximizes the total deviation w.r.t. the 
nominal values, while respecting the bounds on the number of deviations of each 
band. 

We now prove that the robust counterpart of (MILP) under a multi-band sce- 
nario set Si\i can be reformulated as a compact Mixed-Integer Linear Program. 
To this end, we first consider first the linear relaxation (DEVOl-RELAX) of 
(DEVOl) and show that it presents the following nice property. 

Proposition 1. The polytope described by the constraints of (DEVOl-RELAX) 
is integral. 

Proof. We prove the integrality by showing that the coefficient matrix is totally 
unimodular. Consider first the linear relaxation of (DEVOl): 

max E E 4 (DEVOl-RELAX) 

Y^vt^^Ok keK (10) 

E ^ 1 e J (11) 

k£K 

4>0 jeJ,k€K, (12) 
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where we dropped constraints < 1, that are dominated by constraints 
The constraints of (DEVOl-RELAX) have the following matrix form: 



I 



\ 



( : \ 



1 



It is easy to verify that Fi is the incidence matrix of a bipartite graph: the 
elements of the two disjoint set of nodes of the graph are in correspondence with 
the rows of the two distinct layers of block in F,;. Moreover, every column has 
exactly two elements that are not equal to zero: one in the upper layer of blocks 
and one in the lower layer. Being the incidence matrix of a bipartite graph, Fi is 
a totally unimodular matrix |27] . Since Fi is totally unimodular and the vector 
gi is integral, it is well-known that the polytope defined by Fi < gi and jji > 
is integral, thus completing the proof. □ 

Due to this integrality property, we can use strong duality to prove the main 
result of this section. 

Theorem 1. The robust counterpart of a Mixed-Integer Linear Program under 
the multi-band scenario set Sm is equivalent to the following compact Mixed- 
Integer Linear Program: 



max 



E 



5Z '^3 



wf +z^. > dt X 



Ztj > 
Xj > 

Xj e Z4 



'k w. 



'ij 



E 



Zij < bi 



(Rob-MILP) 
i G / 

i e I , j e J, k e K 

i G I,j G J 

i e l,j e K 
je J 

J G Jz C J . 



Proof. As first step, consider the dual problem of (DEVOl-RELAX): 



mm > Uk w, 

keK 



(DEVOl-RELAX-DUAL) 
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Zij > 



j e J,ke K 



where the dual variables wf , Zij are respectively associated with the primal con- 
straints (Uniin]) of (DEVOI-RELAX) defined for constraint i. 

Since (DEVOI-RELAX) is feasible and bounded by definition of the scenario 
set Sm, also (DEVOl-RELAX-DUAL) is feasible and bounded and the optimal 
values of the two problems are the same (strong duality). Then, by Proposi- 
tion [TJ we can replace the maximum total deviation DEVi{x,S]\i) with problem 
(DEVOl-RELAX-DUAL), obtaining the compact Mixed-integer Linear Program 
(Rob-MILP). This concludes the proof. □ 

In comparison to (MILP), this compact formulation uses K ■m + n-m additional 
variables and includes K ■ n ■ m additional constraints. Due to the dominance 
results, we remark that the new formulation (Rob-MILP) has K ■ m variables 
less than the compact counterpart that we presented in [12]. Moreover, we note 
that when (MILP) does not include integer variables (i.e., Jz = 0) and is thus a 
pure Linear Program, (Rob-MILP) is a pure Linear Program as well. 

3 Separation of Robustness Cuts 

As an alternative to directly solving the compact robust counterpart (Rob- 
MILP), we have investigated the development of a cutting-plane algorithm: we 
start by solving the nominal problem (MILP) and then test if the optimal so- 
lution is robust. If not, we generate a cut that imposes robustness [robustness 
cut) and we add it to the problem. This initial step can then be iterated as in a 
typical cutting plane method [57] ■ 

Testing whether a solution x £ M" x zl^'' is robust, corresponds to 
confirming that X^jgj ^ij + DEV[{x,Sm) < h for every constraint i e /. In 
the case of the Bertsimas-Sim model, this is very simple and merely requires to 
sort the deviations and choose the largest Fi values ^T\. In the case of a multi- 
band scenario set, this simple approach docs not guarantee the robustness of a 
computed solution. However, in the following theorem we prove that testing the 
robustness of a solution corresponds to solving a min-cost flow problem [T]. We 
recall that an integral min-cost flow of value n can be computed in polynomial 
time [T]. 

Theorem 2. Let x G M" x zlj;^'' and let Sm &e a multi-band scenario set. 
Moreover, let (G,c)^ be the min-cost flow instance corresponding to x and a 
constraint i £ I of (MILP) and built according to the rules specified in the proof. 
The solution x is robust for constraint i w.r.t. Sm if and only if 



a[x — c*(x) < bi 



where c*{x) is the value of a min cost flow in the instance (G,c)\.. 
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Proof. We prove the statement by showing that there exists a one-to-one cor- 
respondence between the scenarios S G Sm satisfying the profile and bound 
properties (see Lemma [1]) and the integral flows f £ Fi associated with the 
min-cost flow instance (G, c)^ (dcflncd below), with cost relation 

Thus, 

DEV^{x, Sm) = max ^ df^Xj = - min c(/) = ~c*{x), 

which concludes the proof. 

The min-cost flow instance (G, c)^ is deflned as follows. The directed graph 
G contains one vertex Vj for each variable index j S J, one vertex Wk for each 
band k Cz K and two vertices s,t that arc the source and the sink of the flow, 
i.e. V = [Jj^j{vj} U UfeeA'l^'^} ^ i^' 0- The set of arcs A is the union of three 
sets Ai, A2, A3. The set Ai contains one arc from s to every vertex Vj, i.e. 
Ai = {(s, Vj) : j € J}. The set A2 contains one arc from every vertex Vj to every 
vertex Wk, i.e. A2 = {{vj, Wk) '. j (z J,k E K}. Finally, the set A3 contains one arc 
from every vertex Wk to the sink t, i.e. A3 = {{wk^t) : k £ K}. By construction, 
G(V, A) is bipartite and acyclic. To each arc a G ^, we associate an upper bound 
Ua on the flow that can be sent on a and a cost Ca of sending one unit of flow on 
a. These two values are set in the following way: (1,0) when a = {s,Vj) £ Ai; 
{1,—d^jXj) when a = [vj.Wk) G ^2; (^A;,0) when a = {wk,t) G A3. Finally, the 
amount of flow that must be sent trough the network from s to i is equal to n. 
We denote by Fi the set of feasible integral flows of value n of the min-cost flow 
instance (G, c)^. The cost of an (s, t)-flow f £ Fi is deflned by c(/) = J2aeA '^a fa- 
Given a scenario S G Sm that satisfies the profile and bound properties, 
we construct an integral flow / for the min-cost flow instance (G,c)^, in the 
following way: 

= 1 G J (13) 

= 1 ^ 4 = 4 e ^ e ^' (14) 

^0k yke K. (15) 

The flow / respects the arc capacities and has a value of n. Since each deviation 
is assigned by S to exactly one band, / also satisfies fiow conservation in each 
vertex Vj. Furthermore, S assigns exactly 9k values to band fc, and thus 

X! fivj,Wk) = \ Jik{S)\ =0k= f{wk,t) ■ 

Hence, / is in Fi and by definition of the arc costs —c{f) = J2j£j '^fj^j- 

Given a feasible integral flow f G Fi oi value n and cost c(/) of the min-cost 
flow instance (G,c)!p, we construct a scenario S G Sm satisfying the profile and 
bound properties in the following way: 

4 = 4 ^ fiv,,w^) = 1 \/j eJ,ke K. 
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This definition is vahd, since each arc {s,Vj) carries one unit of flow and thus 
due to flow conservation exactly one arc {vj,Wk) carries flow. Furthermore, since 
SfcGif = n-^ exactly 9k arcs with positive flow value enter Wk and therefore 
\JikiS)\ = Ok. Summing up. 

This concludes the proof. □ 

From Theorem [21 we can immediately derive the following corollary: 

Corollary 1. A solution zlj^^' is robust w.r.t. Sm if and only if 

a'^x — c*{x) < bi Vi G /, 

where c*{x) is the value of a minimum cost flow in {G,c)\ . 

According to this corollary, we can test the robustness of a solution x G x 
z|^^' by computing an integral flow of minimum cost c*{x) in (G, c)\, for every 
i G /. If a'x — c*{x) < bi for every i E I, then x is robust feasible. Otherwise x 
is not robust and there exists an index i G / such that a'x — c*{x) > bi and thus 

4 ^3 f{v,.w,) < (16) 

jeJ jeJ keK 

is valid for the polytope of the robust solutions and cuts off the solution x. 



4 Binary Programs with Cost Uncertainties 

In this section, we present results on a subfamily of Mixed-Integer Linear Pro- 
grams, namely pure Binary Programs, and consider the case in which uncertainty 
only affects the objective function. Specifically, we investigate Binary Programs 
of the form: 

min ^CjXj (BP) 

xex c {0,1}" 

with non- negative cost vector, i.e. cj > 0, for ah j G J = {l,...,n}. Rele- 
vant problems like the minimum spanning tree problem, the maximum weighted 
matching problem and the shortest path problem belong to this class of prob- 
lems. 

We study the robust version of (BP) under multi-band uncertainty affecting 
just the cost vector. More formally, for each element j G J we are given the 
nominal cost Cj and a sequence of + 1 deviation values dj, with k G K = 
{0, . . . , A'+}, such that — d^ < d^ < . . . < d^^ < oo (note that in contrast to 
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Section[5]we consider only positive deviations). Through these values, we define: 
1) the zero-deviation band corresponding to the single value = 0; 2) a set 
of positive deviation bands, such that each band k S isr\{0} corresponds to 
the range {d'^^'^,dj]. Furthermore, integer values lk,Uk S Z, < Ik < Uk < n 
represent the lower and upper bounds on the number of deviations falling in 
each band k E K. 

Since (BP) is a special case of (MILP), we can define its robust counterpart 
(Rob-BP) by applying Lemma [1] and Theorem [1] 

min Cj Xj + 9^ Wk + Zj {Rob-BP) 



jeJ keK jej 

Wk + Zj > d*^ Xj j e J,k e K (17) 

Wfc > keK 

Zj >0 j e J 

X G X, 



In order to get a robust optimal solution, we can solve (Rob-BP) by means of 
a commercial MIP solver. However, as an alternative, in what follows we prove 
that we can get a robust optimal solution by solving a sequence of nominal 
problems (BP) with modified objective function. 

Remark 1. Without loss of generality, we can assume that the number of bands 
K is constant and such that K < n. Indeed, if K > n then K' ~ K ~ n 
bands have profile 9k = Q and no coefficient is assigned to them, so they can be 
eliminated from the problem. 

As first step, in the following lemma we prove that for a fixed vector w > 0^, 
solving (Rob-BP) corresponds to solving (BP) with modified cost coefficients. 

Lemma 2. Given w > 0^ , solving (Rob-BP) with w ^ w is equivalent to solving 
the following problem: 

min Cj Xj ( Rob-BP{'w) ) 

X G X , 

with 

Cj = Cj + dj yj G J , 

d, = max jo, max jdf — Wk}} 
■' keK ■' 

Proof. If w = w then (Rob-BP) reduces to: 

J2^kWk+ min Cj Xj + Zj 

Zj > dj Xj — Wk 
z,>0 
X € X, 



j e J,k e K 



(18) 
(19) 
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Because of the binary nature of Xj, constraints (|18m9p can be substituted by 



Zj > dj Xj since: 



Wk, "ik € K 



' " " Zj >niax{0, max{d*^-Wfe}}xj = dj ^j 

We note that dj Xj > by the definition of dj . Moreover, we note that an optimal 
solution (x,z) is such that Zj = djXj, Vj G J. Indeed, suppose by contradiction 
that w.l.o.g. there exists a unique i & J : zi > di xi. Then we can define a 
feasible solution {x' , z') by setting: 

x'j = Xj Vj e J 

4 = ^0 Vj g J\{n 

z[ = Zg — s 

with s = zi, — dgxg > 0, thus obtaining a reduction s in cost and contradicting 
the optimality of (x, w). As Zj = dj Xj at the optimum, for all j G J, we obtain: 

^ ^fc TDfc + min ^ Cj + ^ x^- 
X G X, 

thus completing the proof. □ 

From this lemma it follows immediately that if w* belongs to some optimal 
solution of (Rob-BP), we can solve (Rob-BP(ii;*)) and obtain an optimal robust 
solution X* . 

Corollary 2. If (x* , z* ,w*) is an optimal solution to (Rob-BP) and Z* is the 
corresponding optimal value, then Z* is such that: 

Z* w*k + min ^ Cj Xj . 

keK je.J 

Proof. Let and J^{w) be the sets of feasible solutions of (Rob-BP) and (Rob- 
BP (w)) for some fixed w > 0^ , respectively. Then 

Z* — min > Cj x, + > ^?fc Wfe + > Zj 
^ ■ ' ' jej keK jeJ 

JGJ keK j£j 

keK ^ ' ' ^ ' jeJ jeJ 



^0kwl + unnY^ Cj Xj + ^ dj 

keK jeJ jeJ 

Ok wl + uim^CjXj 
keK ^ jeJ 
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where the equality (a) holds since Zj — djXj at the optimum, as showed in the 
proof of Lemma [21 □ 

Using the previous result, a naive approach to solve (Rob-BP) would be to test 
the value of (Rob-BP(w)) for every w > 0^. However, this would require an 
exponential number of tests. We prove instead that we can limit our attention 
to a subset of vectors w of polynomial size, as indicated in the algorithm and 
theorem that follows. 



Algorithm 1.1 Computing an optimal robust solution of (Rob-BP) 
Input: set of feasible solutions X C {0, 1}"; 
cost values Cj > for all j € J; 

deviation values = dj' < dj < d| < . . . < df for all j G J; 
lower and upper bounds Ik, Uk for all k £ K. 

Output: an optimal robust solution x* £ X 

Step 1: a) set Z — cxj, (x, w, z) — NULL 

b) define Bky = {d^ ~ d']' G J} U {0}, yk, k' € K : k ^ k' 

c) define Ck = {d*^,Vj £ J} U {0}, V/t G K 

Step 2: FOR each possible combination of values b^k' £ B^^' , Cfc G C'k compute a 
solution w to the linear system 

Wk — Wk' ^ bfefc' k,k' G K : k ^ k' 

Wk < Ck k e K 

Wk>Q keK 

IF a solution w exists THEN 

a) solve (Rob-BP(u')) and let {x{w),'W, z{w)), Z(w) be the corresponding 
optimal solution and value, respectively 

b) IF e'w + Z{w) < Z THEN 

i) set Z = O'w + Z(w) and (x, w, z) = (x{'w),w, z(w)) 

Step 3: return {x, w, z) 



Theorem 3. Algorithm 1.1 computes a robust optimal solution to (Rob-BP) 
by solving at most [n -f 1)*^ nominal problems with modified cost coefficients 
(Rob-BP(w)). If there exists an algorithm to solve the nominal problem (BP) in 
polynomial time, Algorithm 1.1 has a polynomial run-time for a constant number 
K of deviation bands. 

Proof. Let x* be an optimal robust solution and let S' be a scenario in Sm such 
that the deviation on x* is maximized. We then define a partition of J according 
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to the bands in which every deviation value d^{j) fahs, i.e. Jk ~ {j € J ■ dj G 
{d^~^, dj]} (this is a partition as every coefficient deviates in exactly one band). 

Claim 1. Let x* be a robust optimal solution. Then there exists (t«*,z*) such 
that (x*, w*, z*) is an optimal solution of (Rob-BP) and w* is a feasible solution 
of the following system LP{bkk', Cfc): 

wl - wl, > bkk' k,k' eK:k^k' {LP{bkk',Ck)) 

Wk < Ck k G K 

wl>0 keK 



with Ck = mmj(=j^{dj X*}, for all k G K and bkk' = max^gj^, {dj' x* - dj x*}, 
for all fc, k' e K : k ^ k' . Moreover, z* = d'j x* - wl, for all k € K,j e Jk- 



To prove the claim, we first note that since x* is a robust optimal solution, 
the corresponding values w* , z* minimize the following linear program obtained 
from (Rob-BP) for x = x*: 

min ^ 9kWk + ^ z^ 
k&K jeJ 

Wk + Zj >djX* j e J,k e K 

Wk>0 k e K 



Zj >0 j G J . 



Its dual problem is: 



max X! X! 4 yj'^ 

^ Vjk <0k kGK 

X! ^1 ^ ^ 

yjk>0 j e J,keK . 

This is exactly the maximum deviation problem DEVOl (sec Section [2]) defined 
for (BP) subject to cost deviation. So an optimal solution y* represents a max- 
imum deviation scenario for x* and we have j G Jk <^=^ y*k = Ij for all k E K. 

Consider now a dual optimal binary solution y* and a primal feasible solution 
(u), z). By the well-known complementary slackness conditions, (it), z) is optimal 
if and only if: 

y*k {wk + Zj - 4 X*) =0 j GJ,keK 

Wk I 0k-y.yL 1=0 k€K 
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Replacing (ui, z) with the inequalities defining any feasible solution {w, z) we 
obtain: 

Vjk (wk + z, - d) x*) = j^J.k^K (20) 

I ^fe-y^yL =0 keK (21) 




= i G J (22) 

Wk + Zj > 4 X* j eJ,k&K (23) 

Wfc > keK (24) 

> j G J (25) 

W.l.o.g. we can assume that the optimal solution y* is such that Ujk ~ 

and J2keK Vjk — ^- j^^t need to consider (in) equalities (PU)) and (23ll25p . 

as the remaining equalities are already satisfied by the choice of y* . Hence, {w, z) 
is optimal if and only if: 

Wk + Zj = dj X* k € K,j e Jk (26) 

Wk + z, > 4 X* keKJe J\Jk (27) 

Wk > k K 

Zj > j & J ■ 

Note that 1^,1^ derive from 1^,^ as y*^ = 1 for all k G K,j G Jk- 
Since the sets Jk, k G K form a partition of J, by (j26p we can set zj = dj x*—Wk 
for all k € K,j G Jk- Yet, we must guarantee that Zj > and so we must have 
Wk < ^ ^ -^^'i S "^fc' which is equivalent to Wfc < miujgjj. dj x*, for 

ah fc G iC. 

Summarizing, for all fc G we have the (in)equalities: 

Wk < min 4 a;* (28) 

= 4 - Wk j G Jfc . (29) 

Consider now an inequality ([27)) and denote by k' the deviation band of coeffi- 
cient j according to x* (i.e., j G Jk')- We can use equality ((29)) in (|27)) obtaining 
■u^fc + (d^' X* - Wk') > d'' X* and thus: 



- Wfc' > dj Xj - Xj 

We can repeat this reasoning for every inequality ()27[) and reorganize the indices 
of the constraints obtaining the system: 

Wk - Wk' > 4 a-'j - df X* k,k' e K : k ^ k' , j G Jk' 
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which is equivalent to 

Wk - Wk' > max{d^' x* - df x*} k,k' e K : k ^ k' 

In total we obtain: given a robust optimal solution x* , (w^z) is robust optimal 
as well if and only if: 

Wk - Wk' > max{d^ x* - df x*} k,k' € K : k k' 

Wk < min d'^ x* fc G K 
jeJk ■' 

Wk>0 k e K 

= X* - Wk k e K,j e Jk ■ 

thus proving Claim 1. 

A crucial consequence of the proved Claim 1 is that we can compute an op- 
timal robust solution x* to (ROB-BP) by solving LP{bkk' ,Ck) for all possible 
combinations of values bkk' , Ck '■ 

bkk' e Bkk' = {4 - 4' , Vj e J} u {0} \^k,k' eK -.k^k' 
Ck e Ck = {d^, Vj e J} u {0} Vfc e K 

Indeed 

1. if w is a feasible solution to LP{bkk' ,Ck), then by Lemma[5]we solve (ROB- 
BP(?Ii)) and obtain the optimal value Z{w) and an optimal solution x{w) of 
(ROB-BP) for fixed w = w; 

2. let wibkk'iCk) be a feasible solution to LP{bkk' ,Ck)- If x(w(bkk' , Ck)) corre- 
sponds with the smallest value Z{w{bkk' ,Ck)) among all the combinations 
i>kk' G Bkk'^Ck € Ck, then x{w{bkk' ,Ck)) is a robust optimal solution to 
(ROB-BP). 

Algorithm 1.1 formally translates these observations by looking for the combina- 
tion bkk' £ Bkk' J Ck G Ck that is associated with the smallest value Z{w{bkk' , c^)). 
Since the number of combinations is equal to {n + 1)'' , we must solve at most 
(n -I- l)'^ problem (ROB-BP(w)) (the problem must not be solved when the 
corresponding LP{-) does not admit a feasible solution w). 

Note that the tractability of the entire algorithm depends on the tractability 
of the nominal optimization problem (BP). If there exists a polynomial algorithm 
for (BP), then Algorithm 1.1 is polynomial as well. 

This completes the proof. □ 

Remark 2. We note that Algorithm 1.1 can be transferred to a robust a- approx- 
imation algorithm if instead of solving (ROB-BP (w)) to optimality we compute 
an a-approximation of the nominal problem. 
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5 Defining a probability bound of constraint violation 



In this section, we analyze the probabihty that a robust optimal solution x* be- 
comes infeasible: a;* is indeed completely protected against deviations captured 
by the multi-band scenario set, but it may still become infeasible under devia- 
tions that lie outside the set. In contrast to the bounds defined in and [8], 
which use a priori information about the not completely known distribution of 
the deviations (specifically, they assume symmetric distributions with symmetric 
range), we define a bound that exploits historical data available on the uncer- 
tain coefiicients. Our basic assumption is that for each coefficient a^ , we have 
W > samples ajj, afj,..., , with each sample lying in the overall deviation 
range, i.e. afj S [ciij — dfj'^^aij + dfj^], a = 1, . . . , W. This is a very reasonable 
assumption in real-world problems, where data about the past behavior of the 
uncertainty arc commonly available. Through this samples, we can define the 
sample mean jiij = '^Ij ^'^"^ define a new typology of proba- 

bilistic bound, as formalized in the proposition that follows. We note that our 
bound holds with a bounded probability, since it is defined using the well-known 
Hoeffding's Inequality 22 : Hoeffding's Inequality bounds the probability that 
the difference between the sample mean and the actual mean of a distribution 
exceeds a fixed threshold. 

Proposition 2. Consider the robust optimization problem (Rob-MILP) under 
a multi-band uncertainty set defined as in Section{l\ and let x* be a robust op- 
timal solution to (Rob-MILP). Suppose that the uncertain coefficients Oij are 
independent random variables and that for each Oij , a number W > Q oj samples 



of, e a, 



with a ~ 1, . . . ,W is available and denote by fiij 

\x* 



the corresponding sample mean. 

Given a constraint i £ I and a constant t > 0, the probability P^{x*) that the 
robust optimal solution x* violates i (i.e., P^ (x*) = P[a^x* > 6,;]) is such that: 



PY{x*) < cxp I In B,j[t,x 



(30) 



where Bij [i, x* 



-a., + df+ - { f,,, + X* (d^+ + rf^-) 



1 



2W\n(3i. 



exp{tx*{a,j -d^j )) - 



-[iM, + x; (4+ + 4-) 



1 



exp {tx*ia,j+df+)) 



Moreover, for fixed r > 0, the bound LSOp holds with probability greater than or 
to 

( 2t^W 
1 — cxp — - 
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Proof. In order to prove the result, we will use the well-known Markov's In- 
equality and Hoeffding's Inequality [35]. Given x* and t > 0, the probability of 
violation of constraint i is: 



P 



aij X* > bi 



(16) 



< cxp(-t&,) Y[E[cxp{ta,jX*)] , (31) 

J 6./ 



where inequality (lb) derives from Markov's Inequality and the independence 
of aij. We are then interested in bounding the moment generating function 
£^ [exp (<aij cc*)] in the r.h.s. of (lb), using the available samples af,-, a = 
1, . . . , W. To this end, we can consider the following formula, that exploits the 
convexity of the exponential function and holds for given 7 > for a bounded 
random variable V : I < V < u [221 : 



£^[e^^] < ^7 [{u-E[V]) eT' + iE[V] - I) e''" 



(32) 



To use this formula, since we do not know the real mean i?[a.y ], in the following 
passages we derive a bound on E[aij]. 

As first step, for each coefficient aij, we use Hoeffding's Inequality to bound 
the probability that the difference between the sample mean and the actual mean 
is above a value r > 0, namely: 



P [x*jE[atj] - Xjfiij >t]=P 



E[aij] ~ ^^^J > — 



< 



By letting (3ij = exp 



< exp 



exp 



2r-'W 



2^ W 



I 



W[x*{a,,+df^ + )~x*{a,, 
2t^W \ 

and by reorganizing the equality through 



simple algebra operations we get t = x*' (d^^ -1- df^ ) \J 2w\ni3 ~ ^^'^ finally 

E[a,,]-li,, > X* {df++df-) 



P 



1 



2Whil3i. 



(33) 



From inequality (1331) , we can derive the following bound for the (unknown) actual 
mean of a,, : 



E[a,j] <^l,,+ X* {df++dt^-) 



1 



2W\nl3i. 



(34) 
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that holds with probability at least 1 — Pij . 

Wc can then use the bound (j34p on the actual mean in the bound (j32p on 
the moment generating function, to finally reach our objective, namely defining 
a bound on E [exp (taij a;*)] . fn our case, the adaptation of (j32p is: 

■ { («u + - E[a^j]) exp {tx*{a,, - dff)) - 

- (a., - df- - E[a,,]) exp {t x* {a,, + df+)) } . (35) 

where we note that we have not yet substituted E[aij\. 

Finally, we substitute the bound ([55]) in (pij) . thus obtaining the bound 
Bij [i, X*] of the statement of the Proposition, that holds with probability greater 
than or equal to rijej(-'^ ~ Pio)- then have: 

exp(-t&i) J]^ [exp (toij a;*)] < exp (-t 6^) J]^ x*] = 

jeJ jeJ 

= exp I hi + In Bij [i, x* 
ending the proof. □ 
6 Conclusions and Future Work 

In this work, we presented new theoretical results abound multi-band uncertainty 
in Robust Optimization. Surprisingly, this natural refinement of the classical 
single band model by Bertsimas and Sim has attracted very little attention and 
we have thus worked on filling the existent knowledge gap. Our ongoing research 
is currently focused on refining the cutting plane method and intensifying the 
computational experiments to other relevant real- world problems, considering 
realistic instances defined in collaboration with our industrial partners. 
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